Figure 4 - huva
Bonaguro et al.
This Rmarkdonwn document can be used to reproduce the panels in figure 4 of the manuscript: huva: human variation in gene expression as surrogate for gene function. To run this script without any problem of dependencies or conflicts in the installation we suggest to use the docker image we provide (see README file).
1 Figure 4 (Supp. Figures S10, S11)
In this figure we expand the huva approach to the entire trascriptome. We first run an huva experiment for all present genes in the 500FG dataset and later explore it, fist manually and later with the use of a co-expression network analysis (CoCena).
1.1 Loading required packages
library(ggplot2)
library(clusterProfiler) # version 3.12.0
library(combinat)
library(ComplexHeatmap)
library(dplyr)
library(ggplot2)
library(graphics)
library(grDevices)
library(grid)
library(igraph)
library(knitr)
library(pcaGoPromoter.Hs.hg19)
library(pheatmap)
library(purrr)
library(RColorBrewer)
library(stringi)
library(tidyr)
library(tidyverse)
library(utils)
library(ggnetwork)
library(intergraph)
library(MCDA)
library(reactome.db)
library(ReactomePA)
library(biomaRt)
library(org.Hs.eg.db) # version 3.8.2
library(Hmisc)
library(gtools)1.2 Costum Functions
source("./source/costum_functions.R")1.3 Loading of the dataset
The calculation of the huva experiment was performed with a legacy version of huva (v 0.1.1), the function used can be found loaded in the envoiroment (calculate_network). Here we will simply load the calculated dataset to keep calculation time short. The result of the caclulation with the current version of huva is identical.
load("./data/network_data.RData")1.4 Manual analysis
1.4.1 Gene involved in Monocytes biology
We where fist interested in the analysis of genes strongly influencing the number of circulating monocytes. We visualized them as volcano plot according to the result of the huva experiment
1.4.1.1 Plot
vulcano_cell_type(x = "Monocytes (CD14+)", y = "pval Monocytes (CD14+)") +
ylab("-log10(p value)") +
theme(legend.position = "none") +
ggtitle("Fig. 4b - Volcano Monocytes")
part1 <- list()
part1$CD14_spec <- intersect(row.names(FG500_bin$cell_fc[FG500_bin$cell_fc$`Monocytes (CD14+)`<1,]),
row.names(FG500_bin$cell_pval[FG500_bin$cell_pval$`Monocytes (CD14+)`<0.05,]))
top_genes_plot(data =log2(FG500_bin$cell_fc)*-log10(FG500_bin$cell_pval),
list = part1$CD14_spec,
paramiter = "Monocytes (CD14+)",
fill = "Monocytes (CD14+)",
title = "Fig 4c - Top 20 - Monocytes",
top_n = 20, decreasing = F, log_trans=F) +
theme(legend.position = "none") +
ylab("log2(FC)*-log10(p val)") +
xlab("")
1.4.2 Gene involved in IFNg secretion
1.4.2.1 Plot
vulcano_cytokine_type(x = "IFNy_PHA_WB_48h", y = "pval IFNy_PHA_WB_48h") +
ylab("-log10(p value)") +
theme(legend.position = "none") +
ggtitle("Fig. S10b - Volcano IFNg")
part1$IFN_spec <- intersect(row.names(FG500_bin$cyto_fc[FG500_bin$cyto_fc$IFNy_PHA_WB_48h<1,]),
row.names(FG500_bin$cyto_pval[FG500_bin$cyto_pval$IFNy_PHA_WB_48h<0.05,]))
top_genes_plot(data =log2(FG500_bin$cyto_fc)*-log10(FG500_bin$cyto_pval),
list = part1$IFN_spec,
paramiter = "IFNy_PHA_WB_48h",
fill = "IFNy_PHA_WB_48h",
title = "Fig. S10c - Top20 - IFNg",
top_n = 20, decreasing = F, log_trans=F) +
theme(legend.position = "none") +
ylab("log2(FC)*-log10(p val)") +
xlab("")
1.5 CoCena Cell counts
We will now explore the result of the transcriptome-wide huva expreiment we used the Co-expression network analysis tool CoCena2. We generate first a network of huva experiments for the cell counts to identify main regulators for each cell type.
1.5.1 Data Preparation
We here load the previously calculated huva experiment results. As robust metric for the calculation of the network we used the product of the log2 fold change and the negative log 10 p value fo the comparison.We also inverted the sign of the result, as described in the manuscript, to ease the interpretation of the results (a negative values represent a negative impact on a specific cell type/cytokine).
load("./data/network_data.RData")
df <- -(log2(FG500_bin$cell_fc)*-log10(FG500_bin$cell_pval))
count_file_name <- df
count_file_name <- count_file_name[, c(1,14,27,33,53,59)]
info_dataset_name <- data.frame(gene=colnames(count_file_name))
info_dataset_name$data_type <- "cell_count"
rownames(info_dataset_name) <- info_dataset_name$gene
info_dataset_name$group <- info_dataset_name$gene
rm(df)1.5.2 Settings of the analysis
working_directory = "./"
topvar_genes = 3000
voi = "group"
TF_list_name = "TFcat.txt"
gmtfile_name_hallmarks = "h.all.v6.1.symbols.gmt"
gmtfile_name_go = "c5.bp.v7.0.symbols.gmt"
gmtfile_name_kegg = "c2.cp.kegg.v7.0.symbols.gmt"
gmtfile_name_reactome = "c2.cp.reactome.v7.0.symbols.gmt"
organism = "human"
min_corr=0.5
range_cutoff_length=1000
print_distribution_plots = FALSE
min_nodes_number_for_network=20
min_nodes_number_for_cluster=20
data_in_log=T
range_GFC=2.0
layout_algorithm = "layout_with_fr" 1.5.2.1 References for functional enrichment
mart <- biomaRt::useMart("ensembl")
mart <- biomaRt::listDatasets(mart)
human <- biomaRt::useMart("ensembl", dataset = "hsapiens_gene_ensembl")
mouse <- biomaRt::useMart("ensembl", dataset = "mmusculus_gene_ensembl")
count_table <- count_file_name
universe_Entrez <- clusterProfiler:: bitr(row.names(count_table),
fromType="SYMBOL",
toType="ENTREZID",
OrgDb="org.Hs.eg.db",
drop = T)
info_dataset <- info_dataset_name
TF_list <- read.delim(paste0(working_directory, "reference_files/", TF_list_name),
header=TRUE,
check.names=F)
gmtfile_hallmarks <- clusterProfiler::read.gmt(paste0(working_directory, "reference_files/", gmtfile_name_hallmarks))
gmtfile_go <- clusterProfiler::read.gmt(paste0(working_directory, "reference_files/", gmtfile_name_go))
gmtfile_kegg <- clusterProfiler::read.gmt(paste0(working_directory, "reference_files/", gmtfile_name_kegg))
gmtfile_reactome <- clusterProfiler::read.gmt(paste0(working_directory, "reference_files/", gmtfile_name_reactome))1.5.3 Filtering the data
We filter tha data here for the top 3000 most variable huva experiments
ds = count_table[order(apply(count_table,1,var), decreasing=T),]
dd2 <- head(ds,topvar_genes)
dd2 = t(dd2)1.5.4 Correlation and cut-off calculation
source(paste0(working_directory,"source/", "correlation_actions.R"))
source(paste0(working_directory,"source/", "obtain_cutoff_stats.R"))
cutoff_stats = do.call("rbind", lapply(X = range_cutoff,
FUN = cutoff_prep,
corrdf_r = correlation_df_filt,
print.all.plots = print_distribution_plots))
source(paste0(working_directory,"source/", "optimal_cutoff.R"))## [1] "The calculated optimal cutoff is: 0.971"1.5.5 Data filtering based on correlation cut off
optimal_cutoff = calculated_optimal_cutoff
cutoff_wd <- paste0("Cell_count_",optimal_cutoff, "_", topvar_genes)
if(!cutoff_wd %in% list.dirs(working_directory)) {
dir.create(paste0(working_directory,cutoff_wd))}## Warning in dir.create(paste0(working_directory, cutoff_wd)): './
## Cell_count_0.971_3000' already existsstats_optimal_cutoff <- cutoff_stats[cutoff_stats$cutoff == optimal_cutoff, c("degree", "Probs")]
filt_cutoff_data = correlation_df_filt %>% dplyr::filter(rval > optimal_cutoff)
filt_cutoff_graph = igraph::graph_from_data_frame(filt_cutoff_data,directed=FALSE)
filt_cutoff_counts = ds[row.names(ds) %in% names(V(filt_cutoff_graph)),]
corresp_info = info_dataset[rownames(dd2)%in%rownames(info_dataset),]
print(paste("After using the optimal cutoff of",optimal_cutoff, "the number of edges =",
nrow(filt_cutoff_data), "and the number of nodes =", nrow(filt_cutoff_counts)))## [1] "After using the optimal cutoff of 0.971 the number of edges = 61585 and the number of nodes = 2894"1.5.6 GFC calculation
source(paste0(working_directory,"source/", "GFC_calculation.R" ))
GFC <- GFC_calculation(voi_id = voi)
GFC_all_samples <- GFC_calculation(voi_id = "gene")1.5.7 Fig 4d - Clustering
source(paste0(working_directory,"source/", "cluster_calculation.R" ))
cluster_information <- cluster_calculation(igraph = filt_cutoff_graph,
cluster_algo = "auto",
no_of_iterations = 10,
max_cluster_count_per_gene = 10,
min_cluster_size = min_nodes_number_for_cluster,
GFC = GFC)
cluster_information_all <- cluster_calculation(igraph = filt_cutoff_graph,
cluster_algo = "auto",
no_of_iterations = 10,
max_cluster_count_per_gene = 10,
min_cluster_size = min_nodes_number_for_cluster,
GFC = GFC_all_samples)source(paste0(working_directory,"source/", "heatmap_clusters.R" ))
heatmap_cluster <- heatmap_clusters(data = cluster_information, cluster_cols = T)
1.5.8 Network generation
source(paste0(working_directory,"source/", "network_generation.R" ))
return_network <- network_layout(igraph.object = filt_cutoff_graph,
use.layout = layout_algorithm,
min.nodes.number = min_nodes_number_for_network)
igraph_object <- return_network$graph_object
layout_matrix <- return_network$layout
node_attributes <- node_information(igraph.object = igraph_object,
data_df = cluster_information,
GFC_df = GFC,
TF_df = TF_list,
hallmark_df = gmtfile_hallmarks,
go_df = gmtfile_go,
kegg_df = gmtfile_kegg,
reactome_df = gmtfile_reactome,
org = organism)
network_object <- generate_network_object(graph_object = igraph_object,
use_layout = layout_matrix)1.5.9 Fig. 4e/f - Network visualization
source(paste0(working_directory,"source/", "network_visualization.R" ))
network_cluster <- visualize_network(network = network_object,
color.by = "cluster",
select.cluster = NULL,
plot.subnetwork = NULL,
gene.label = NULL,
use.layout = layout_algorithm,
save.pdf=F)
print(network_cluster)
source(paste0(working_directory,"source/", "network_visualization.R" ))
network_GFC <- GFC_colored_network(network = network_object,
select.cluster = NULL,
plot.subnetwork = NULL,
gene.label = NULL,
use.layout = layout_algorithm,
save.pdf = F,
save.single.pdf = F)
gridExtra::marrangeGrob(grobs = network_GFC, ncol = 1, nrow = 1)
1.5.9.1 Fig S12b - STAT1 in the network
source(paste0(working_directory,"source/", "network_visualization.R" ))
STAT1 <- visualize_network(network = network_object,
color.by = "cluster",
select.cluster = c("khaki"),
gene.label = c("STAT1"),
plot.subnetwork = NULL,
use.layout = layout_pers,
save.pdf = F)
print(STAT1)
1.5.10 Fig. 4g - Cluster Profiler
source(paste0(working_directory,"source/","clusterprofiler_autoCena.R"))
clust_prof= clusterprofiler_autoCena(cluster_data = cluster_information,
cutoff_wd = cutoff_wd,
originalwd = working_directory,
chosen_cutoff = optimal_cutoff,
group = voi)
names(clust_prof) <- cluster_information[cluster_information$cluster_included=="yes",]$colorsource("./source/plot_selected_GOEA.R")
1.6 CoCena Cytokines
1.6.1 Data Preparation
load("./data/network_data.RData")
df <- -(log2(FG500_bin$cyto_fc)*-log10(FG500_bin$cyto_pval))
count_file_name <- df
info_dataset_name <- data.frame(gene=colnames(count_file_name))
info_dataset_name$data_type <- "cytokines"
rownames(info_dataset_name) <- info_dataset_name$gene
info_dataset_name$group <- unlist(lapply(strsplit(rownames(info_dataset_name), "_"), `[[`, 1))
info_dataset_name$stim <- unlist(lapply(strsplit(rownames(info_dataset_name), "_"), `[[`, 2))
info_dataset_name$time <- unlist(lapply(strsplit(rownames(info_dataset_name), "_"), `[[`, 4))
rm(df)1.6.2 Settings of the analysis
working_directory = "./"
topvar_genes = 3000
voi = "group"
TF_list_name = "TFcat.txt"
gmtfile_name_hallmarks = "h.all.v6.1.symbols.gmt"
gmtfile_name_go = "c5.bp.v7.0.symbols.gmt"
gmtfile_name_kegg = "c2.cp.kegg.v7.0.symbols.gmt"
gmtfile_name_reactome = "c2.cp.reactome.v7.0.symbols.gmt"
organism = "human"
min_corr=0.5
range_cutoff_length=1000
print_distribution_plots = FALSE
min_nodes_number_for_network=20
min_nodes_number_for_cluster=20
data_in_log=T
range_GFC=2.0
layout_algorithm = "layout_with_fr" 1.6.2.1 References for functional enrichment
mart <- biomaRt::useMart("ensembl")
mart <- biomaRt::listDatasets(mart)
human <- biomaRt::useMart("ensembl", dataset = "hsapiens_gene_ensembl")
mouse <- biomaRt::useMart("ensembl", dataset = "mmusculus_gene_ensembl")
count_table <- count_file_name
universe_Entrez <- clusterProfiler:: bitr(row.names(count_table),
fromType="SYMBOL",
toType="ENTREZID",
OrgDb="org.Hs.eg.db",
drop = T)
info_dataset <- info_dataset_name
TF_list <- read.delim(paste0(working_directory, "reference_files/", TF_list_name),
header=TRUE,
check.names=F)
gmtfile_hallmarks <- clusterProfiler::read.gmt(paste0(working_directory, "reference_files/", gmtfile_name_hallmarks))
gmtfile_go <- clusterProfiler::read.gmt(paste0(working_directory, "reference_files/", gmtfile_name_go))
gmtfile_kegg <- clusterProfiler::read.gmt(paste0(working_directory, "reference_files/", gmtfile_name_kegg))
gmtfile_reactome <- clusterProfiler::read.gmt(paste0(working_directory, "reference_files/", gmtfile_name_reactome))1.6.3 Filtering the data
Also here we take the top 3000 huva experiments
ds = count_table[order(apply(count_table,1,var), decreasing=T),]
dd2 <- head(ds,topvar_genes)
dd2 = t(dd2)1.6.4 Correlation and cut-off calculation
source(paste0(working_directory,"source/", "correlation_actions.R"))
source(paste0(working_directory,"source/", "obtain_cutoff_stats.R"))
cutoff_stats = do.call("rbind", lapply(X = range_cutoff,
FUN = cutoff_prep,
corrdf_r = correlation_df_filt,
print.all.plots = print_distribution_plots))
source(paste0(working_directory,"source/", "optimal_cutoff.R"))## [1] "The calculated optimal cutoff is: 0.674"1.6.5 Data filtering based on correlation cut off
optimal_cutoff = calculated_optimal_cutoff
cutoff_wd <- paste0("Cytokines_",optimal_cutoff, "_", topvar_genes)
if(!cutoff_wd %in% list.dirs(working_directory)) {
dir.create(paste0(working_directory,cutoff_wd))}## Warning in dir.create(paste0(working_directory, cutoff_wd)): './
## Cytokines_0.674_3000' already existsstats_optimal_cutoff <- cutoff_stats[cutoff_stats$cutoff == optimal_cutoff, c("degree", "Probs")]
filt_cutoff_data = correlation_df_filt %>% dplyr::filter(rval > optimal_cutoff)
filt_cutoff_graph = igraph::graph_from_data_frame(filt_cutoff_data,directed=FALSE)
filt_cutoff_counts = ds[row.names(ds) %in% names(V(filt_cutoff_graph)),]
corresp_info = info_dataset[rownames(dd2)%in%rownames(info_dataset),]
print(paste("After using the optimal cutoff of",optimal_cutoff, "the number of edges =",
nrow(filt_cutoff_data), "and the number of nodes =", nrow(filt_cutoff_counts)))## [1] "After using the optimal cutoff of 0.674 the number of edges = 159827 and the number of nodes = 2943"1.6.6 GFC calculation
source(paste0(working_directory,"source/", "GFC_calculation.R" ))
GFC <- GFC_calculation(voi_id = voi)
GFC_all_samples <- GFC_calculation(voi_id = "stim")1.6.7 Clustering
source(paste0(working_directory,"source/", "cluster_calculation.R" ))
cluster_information <- cluster_calculation(igraph = filt_cutoff_graph,
cluster_algo = "auto",
no_of_iterations = 10,
max_cluster_count_per_gene = 10,
min_cluster_size = min_nodes_number_for_cluster,
GFC = GFC)1.6.7.1 Fig. S10d
source(paste0(working_directory,"source/", "heatmap_clusters.R" ))
heatmap_cluster <- heatmap_clusters(data = cluster_information, cluster_cols = T)
cluster_information_all <- cluster_calculation(igraph = filt_cutoff_graph,
cluster_algo = "auto",
no_of_iterations = 10,
max_cluster_count_per_gene = 10,
min_cluster_size = min_nodes_number_for_cluster,
GFC = GFC_all_samples)1.6.7.2 Fig. S11a
source(paste0(working_directory,"source/", "heatmap_clusters.R" ))
heatmap_cluster_all_samples <- heatmap_clusters(data = cluster_information_all, cluster_cols = T)
1.6.8 Network generation
source(paste0(working_directory,"source/", "network_generation.R" ))
return_network <- network_layout(igraph.object = filt_cutoff_graph,
use.layout = layout_algorithm,
min.nodes.number = min_nodes_number_for_network)
igraph_object <- return_network$graph_object
layout_matrix <- return_network$layout
node_attributes <- node_information(igraph.object = igraph_object,
data_df = cluster_information,
GFC_df = GFC,
TF_df = TF_list,
hallmark_df = gmtfile_hallmarks,
go_df = gmtfile_go,
kegg_df = gmtfile_kegg,
reactome_df = gmtfile_reactome,
org = organism)
network_object <- generate_network_object(graph_object = igraph_object,
use_layout = layout_matrix)
node_attributes <- node_information(igraph.object = igraph_object,
data_df = cluster_information,
GFC_df = GFC_all_samples,
TF_df = TF_list,
hallmark_df = gmtfile_hallmarks,
go_df = gmtfile_go,
kegg_df = gmtfile_kegg,
reactome_df = gmtfile_reactome,
org = organism)
network_object_all_samples <- generate_network_object(graph_object = igraph_object,
use_layout = layout_matrix)1.6.9 Fig. S10e/f - Network visualization
source(paste0(working_directory,"source/", "network_visualization.R" ))
network_cluster <- visualize_network(network = network_object,
color.by = "cluster",
select.cluster = NULL,
plot.subnetwork = NULL,
gene.label = NULL,
use.layout = layout_algorithm,
save.pdf=F)
print(network_cluster)
source(paste0(working_directory,"source/", "network_visualization.R" ))
network_GFC <- GFC_colored_network(network = network_object,
select.cluster = NULL,
plot.subnetwork = NULL,
gene.label = NULL,
use.layout = layout_algorithm,
save.pdf = F,
save.single.pdf = F)
gridExtra::marrangeGrob(grobs = network_GFC, ncol = 1, nrow = 1)
1.6.10 Fig. S10g - Cluster Profiler
source(paste0(working_directory,"source/","clusterprofiler_autoCena.R" ))
clust_prof= clusterprofiler_autoCena(cluster_data = cluster_information,
cutoff_wd = cutoff_wd,
originalwd = working_directory,
chosen_cutoff = optimal_cutoff,
group = voi)
names(clust_prof) <- cluster_information[cluster_information$cluster_included=="yes",]$colorsource("./source/plot_selected_GOEA_cyto.R")
2 Clean the environment
source("./source/clean.R")3 Session info
info <- sessionInfo()
info## R version 4.0.3 (2020-10-10)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 20.04 LTS
##
## Matrix products: default
## BLAS/LAPACK: /usr/lib/x86_64-linux-gnu/openblas-pthread/libopenblasp-r0.3.8.so
##
## locale:
## [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
## [3] LC_TIME=en_US.UTF-8 LC_COLLATE=en_US.UTF-8
## [5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=C
## [7] LC_PAPER=en_US.UTF-8 LC_NAME=C
## [9] LC_ADDRESS=C LC_TELEPHONE=C
## [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
##
## attached base packages:
## [1] parallel stats4 grid stats graphics grDevices utils
## [8] datasets methods base
##
## other attached packages:
## [1] GSEABase_1.52.1 graph_1.68.0
## [3] annotate_1.68.0 XML_3.99-0.5
## [5] gtools_3.8.2 Hmisc_4.4-2
## [7] Formula_1.2-4 survival_3.2-7
## [9] lattice_0.20-41 org.Hs.eg.db_3.8.2
## [11] biomaRt_2.46.3 ReactomePA_1.34.0
## [13] reactome.db_1.74.0 AnnotationDbi_1.52.0
## [15] IRanges_2.24.1 S4Vectors_0.28.1
## [17] Biobase_2.50.0 BiocGenerics_0.36.0
## [19] MCDA_0.0.20 intergraph_2.0-2
## [21] ggnetwork_0.5.8 forcats_0.5.1
## [23] stringr_1.4.0 readr_1.4.0
## [25] tibble_3.0.6 tidyverse_1.3.0
## [27] tidyr_1.1.2 stringi_1.5.3
## [29] RColorBrewer_1.1-2 purrr_0.3.4
## [31] pheatmap_1.0.12 pcaGoPromoter.Hs.hg19_1.26.0
## [33] knitr_1.31 igraph_1.2.6
## [35] dplyr_1.0.4 ComplexHeatmap_2.6.2
## [37] combinat_0.0-8 clusterProfiler_3.12.0
## [39] ggplot2_3.3.3
##
## loaded via a namespace (and not attached):
## [1] readxl_1.3.1 shadowtext_0.0.7 backports_1.2.1
## [4] circlize_0.4.12 fastmatch_1.1-0 BiocFileCache_1.14.0
## [7] plyr_1.8.6 splines_4.0.3 BiocParallel_1.24.1
## [10] digest_0.6.27 htmltools_0.5.1.1 GOSemSim_2.16.1
## [13] viridis_0.5.1 GO.db_3.12.1 magrittr_2.0.1
## [16] checkmate_2.0.0 memoise_2.0.0 cluster_2.1.1
## [19] openxlsx_4.2.3 graphlayouts_0.7.1 modelr_0.1.8
## [22] matrixStats_0.58.0 askpass_1.1 enrichplot_1.10.2
## [25] prettyunits_1.1.1 jpeg_0.1-8.1 colorspace_2.0-0
## [28] blob_1.2.1 rvest_0.3.6 rappdirs_0.3.3
## [31] ggrepel_0.9.1 haven_2.3.1 xfun_0.21
## [34] crayon_1.4.1 jsonlite_1.7.2 scatterpie_0.1.5
## [37] Rglpk_0.6-4 glue_1.4.2 polyclip_1.10-0
## [40] gtable_0.3.0 GetoptLong_1.0.5 graphite_1.36.0
## [43] shape_1.4.5 scales_1.1.1 DOSE_3.16.0
## [46] DBI_1.1.1 Rcpp_1.0.6 xtable_1.8-4
## [49] htmlTable_2.1.0 viridisLite_0.3.0 progress_1.2.2
## [52] clue_0.3-58 foreign_0.8-81 bit_4.0.4
## [55] htmlwidgets_1.5.3 httr_1.4.2 fgsea_1.12.0
## [58] ellipsis_0.3.1 pkgconfig_2.0.3 farver_2.0.3
## [61] nnet_7.3-15 dbplyr_2.1.0 labeling_0.4.2
## [64] tidyselect_1.1.0 rlang_0.4.10 reshape2_1.4.4
## [67] munsell_0.5.0 cellranger_1.1.0 tools_4.0.3
## [70] cachem_1.0.4 cli_2.3.0 generics_0.1.0
## [73] RSQLite_2.2.3 broom_0.7.4 evaluate_0.14
## [76] fastmap_1.1.0 yaml_2.2.1 bit64_4.0.5
## [79] fs_1.5.0 tidygraph_1.2.0 zip_2.1.1
## [82] ggraph_2.0.4 slam_0.1-48 DO.db_2.9
## [85] xml2_1.3.2 compiler_4.0.3 rstudioapi_0.13
## [88] curl_4.3 png_0.1-7 reprex_1.0.0
## [91] tweenr_1.0.1 highr_0.8 Matrix_1.3-2
## [94] vctrs_0.3.6 pillar_1.4.7 lifecycle_1.0.0
## [97] BiocManager_1.30.12 GlobalOptions_0.1.2 data.table_1.13.6
## [100] cowplot_1.1.1 qvalue_2.22.0 latticeExtra_0.6-29
## [103] R6_2.5.0 network_1.16.1 gridExtra_2.3
## [106] MASS_7.3-53.1 assertthat_0.2.1 openssl_1.4.3
## [109] rjson_0.2.20 withr_2.4.1 hms_1.0.0
## [112] rpart_4.1-15 glpkAPI_1.3.2 rmarkdown_2.6
## [115] rvcheck_0.1.8 Cairo_1.5-12.2 ggforce_0.3.2
## [118] base64enc_0.1-3 lubridate_1.7.9.24 Save environment
save.image(paste("./data/",Sys.Date(), "_huva_Figure_4.RData", sep = ""))